Climatologically Significant Effects of Some Approximations in the Bulk Parameterizations of Turbulent Air – Sea Fluxes

This paper quantifies the impacts of approximations and assumptions in the parameterization of bulk formulas on the exchange of momentum, heat, and freshwater computed between the ocean and atmosphere. An ensemble of sensitivity experiments is examined. Climatologies of wind stress, turbulent heat flux, and evaporation for the period 1982–2014 are computed using SST and surface meteorological state variables from ERA-Interim. Each experiment differs from the defined control experiment in only one aspect of the parameterization of the bulk formulas. The wind stress is most sensitive to the closure used to relate the neutral drag coefficient to the wind speed in the bulk algorithm, which mainly involves the value of the Charnock parameter. The disagreement between the state-of-the-art algorithms examined is typically on the order of 10%. The largest uncertainties in turbulent heat flux and evaporation are also related to the choice of the algorithm (typically 15%) but also emerge in experiments examining approximations related to the surface temperature and saturation humidity. Thus, approximations for the skin temperature and the saltrelated reduction of saturation humidity have a substantial impact on the heat flux and evaporation (typically 10%). Approximations such as the use of a fixed air density, sea level pressure, or simplified formula for the saturation humidity lead to errors no larger than 4% when tested individually. The impacts of these approximations combine linearly when implemented together, yielding errors up to 20%overmidand subpolar latitudes.


Introduction
Accurate estimation of the exchanges of momentum, heat, and freshwater between the ocean and atmosphere is critical for a wide range of weather-and climate-related studies.The turbulent air-sea fluxes (TASFs), which include wind stress, evaporation, and the latent and sensible heat flux components, are the primary mechanism by which these exchanges occur.As such, TASFs influence the variability and climatology of the net surface heat flux at all scales.In contrast to the radiative shortwave and longwave heat flux components, the spatial pattern of TASFs is characterized by strongly localized maxima over the mid-and high latitudes, making them crucial drivers of the variability of many ocean processes, such as deep convection in the subpolar waters (Visbeck et al. 2003;Moore et al. 2014;Holdsworth and Myers 2015).
TASFs also have a strong impact on the atmospheric circulation and climate.For instance, in the western boundary current regions, TASFs can modulate the atmospheric responses at mesoscales and synoptic scales (Zolina and Gulev 2003;Small et al. 2008;Minobe et al. 2008;Ma et al. 2015b,a).
Since direct measurements of TASFs tend to be idealized, infrequent, and highly localized, they cannot be used to build global or regional climatologies of TASFs (Brunke et al. 2003(Brunke et al. , 2011)).Rather, these data provide valuable constraints for the development and improvement of parameterizations (Fairall et al. 1996b(Fairall et al. , 2003(Fairall et al. , 2011;;Edson et al. 2013).
Hence, global estimates of TASFs are derived using bulk formulas, which relate each TASF to more easily measurable and widely available meteorological surface state variables (SSVs), and a bulk transfer coefficient (BTC).BTCs are determined with what we will refer to as a bulk algorithm.The core of a typical bulk algorithm includes dependencies of BTCs on the wind speed and the stability of the atmospheric surface layer (ASL) as well as the adjustment of atmospheric scalars to the standard height through the flux-profile relationships.Traditionally, bulk algorithms are developed for in situ point measurements of SSVs, with neutral drag, heat, and moisture BTCs derived from ship and buoy measurements (Large andPond 1981, 1982;Smith 1988;Fairall et al. 1996bFairall et al. , 2003;;Bradley and Fairall 2007).In this sense, bulk parameterizations are most suited to applications utilizing voluntary observing ship (VOS) and buoy data only.Nevertheless, for most purposes, bulk formulas are used to compute TASFs using a variety of SSV datasets.
Discrepancies between state-of-the-art bulk algorithms have been repeatedly shown to lead to a relatively large spread in the computed TASFs (Blanc 1985(Blanc , 1986(Blanc , 1987;;Zeng et al. 1998;Eymard et al. 1999;Brunke et al. 2002Brunke et al. , 2003)).Moreover, given the large uncertainties prevalent in SSVs from different sources (Josey et al. 1999;Berry and Kent 2005;Josey et al. 2014), it is common practice to neglect the accuracy of various constants and parameter approximations in bulk formulas.Yet, the uncertainties associated with these assumptions remain too poorly quantified to constitute actionable information for concerned users.An important question is whether or not the uncertainties associated with assumptions and approximations in the bulk formulas are significant for climatological estimates of TASFs.In this study, we address these issues and quantify the effect of various assumptions and approximations in the bulk formulas on the global TASF estimates.This is assessed through numerical experiments in which global TASFs are computed over a set of configurable bulk routines using prescribed daily SST and 3-hourly atmospheric SSVs from the ERA-Interim reanalysis for the period 1982-2014.
By computing TASFs with a fixed bulk formulas frame, modifiable constants, interchangeable bulk algorithms, and various parameter approximations, we quantify errors with respect to a reference control experiment.
This paper addresses the sensitivity of computed TASFs to various aspects of the parameterization of the bulk formulas with respect to the uncertainties arising due to the choice of the bulk algorithm.These aspects include the air density, sea surface saturation humidity, mean sea level pressure, cool skin and warm layer, and sea surface currents.
This paper is organized as follows: In section 2, we provide a short overview of the bulk formulas and some well-established bulk algorithms and describe the setup of the control experiment.In section 3, we quantify TASF uncertainties related to various assumptions and approximations in the bulk formulas, each tested with an individual sensitivity experiment.In section 4, we focus on the uncertainties related to the choice and parameterization of the bulk algorithm used to compute the BTCs.Our conclusions are summarized in section 5.
A detailed nomenclature relating the acronyms and symbols used in this paper to the terms or physical parameters they stand for is outlined in appendix A.

Method and experiment design
a.The bulk model Turbulent air-sea fluxes are computed using the sea surface properties and atmospheric SSVs at height z above the sea surface, with the traditional aerodynamic bulk formulas: with u z ' T z 1 gz, and (1e) where t is the wind stress, Q H is the sensible heat flux, E is the evaporation, and Q L is the latent heat flux.
Throughout this paper, we use the convention that a positive sign of t, Q H , and Q L means a gain of the relevant quantity for the ocean, while a positive E implies a freshwater loss for the ocean.The term r is the density of air; C D , Q H , and C E are the BTCs for momentum, sensible heat, and moisture, respectively; C P is the heat capacity of moist air, and L y is the latent heat of vaporization of water; u z , T z , and q z are the potential temperature, temperature, and specific humidity of air at height z, respectively; gz is a temperature correction term, which accounts for the adiabatic lapse rate and approximates the potential temperature at height z (Josey et al. 2013); and U z is the wind speed vector at height z (possibly referenced to the surface current u 0 ; section 3e).The bulk scalar wind speed U B is the scalar wind speed jU z j with the potential inclusion of a gustiness contribution (section 4c).The term P 0 is the mean sea level pressure (SLP), T s is the sea surface temperature, and q s is the saturation-specific humidity of air at temperature T s and includes a 2% reduction to account for the presence of salt in seawater (Sverdrup et al. 1942;Kraus and Businger 1996).Depending on the bulk parameterization used, T s can be the temperature at the airsea interface [sea surface skin temperature (SSST)] or at a few tens of centimeters below the surface [bulk sea surface temperature (SST)].The SSST differs from the SST due to the contributions of two effects of opposite sign: the cool skin and warm layer (CSWL).The cool skin refers to the cooling of the millimeter-scale uppermost layer of the ocean, in which the net upward flux of heat to the atmosphere is ineffectively sustained by molecular diffusion.As such, a steep vertical gradient of temperature must exist to ensure the heat flux continuity with underlying layers in which the same flux is sustained by turbulence.The warm layer refers to the warming of the upper few meters of the ocean under sunny conditions.The CSWL effects are most significant under weak wind conditions due to the absence of substantial surface vertical mixing (caused by, e.g., breaking waves).The impact of the CSWL on the computed TASFs is discussed in section 3d.Accuracy of the bulk formulas strongly relies on the BTCs C D , C H , and C E .In a typical bulk algorithm, the BTCs under neutral stability conditions are defined using in situ flux measurements, while their dependence on the stability is accounted through the Monin-Obukhov similarity theory and the flux-profile relationships (e.g., Paulson 1970).BTCs are functions of the wind speed and the near-surface stability of the ASL and hence depend on U B , T s , T z , q s , and q z .
In this study, we focus on three of the most common algorithms used by the GCM community: COARE, after Fairall et al. (1996bFairall et al. ( , 1997Fairall et al. ( , 2003Fairall et al. ( , 2011) ) and more recently Edson et al. (2013); NCAR, following the work of (Large andPond 1981, 1982;Large et al. 1997;Large andYeager 2004, 2009); and ECMWF (Miller et al. 1992;Beljaars 1995Beljaars , 1997;;Zeng and Beljaars 2005).By COARE, we explicitly refer to the COARE 3.0 algorithm (Fairall et al. 2003), by NCAR we refer to the algorithm developed by Large and Yeager (2009), and by ECMWF we refer to the version of the bulk algorithm used in recent cycles of the Integrated Forecast System (IFS) of the ECMWF, such as cycle 40 (ECMWF 2014).

b. Experimental setup
We use bulk formulas with 3-hourly atmospheric SSVs and daily mean SST from ERA-Interim (Dee et al. 2011) to compute global 3-hourly estimates of the TASFs from 1982 to 2014.Our experiments are initiated in 1982 as the SST data product used as the surface boundary condition in ERA-Interim becomes truly daily from the end of 1981.Prior to this, daily mean SSTs were derived through interpolation of monthly mean data.
In the control experiment CtrlE (Table 1), TASFs are computed at the standard height of z 5 10 m, according to Eqs. (1) and the procedure outlined in Fig. 1.BTCs are computed using the COARE 3.0 algorithm, which has consistently performed well in algorithm intercomparisons (e.g., Brunke et al. 2003).Sensitivity experiments using the NCAR and ECMWF bulk algorithms are also performed to assess the algorithm-related TASFs uncertainties (section 4).
The use of the COARE algorithm in CtrlE implies that T s is the SSST.Therefore, the SST is corrected online to account for the CSWL effects (Fairall et al. 1996a).The CSWL correction requires knowledge of the solar and the nonsolar components of the surface net heat flux to the ocean: with where a and d are the albedo and emissivity of the sea surface, respectively.Thus, we use the computed Q L and Q H and the 3-hourly surface downwelling shortwave and  longwave radiative fluxes (Q swY and Q lwY , respectively) from ERA-Interim to correct the daily SST every 3 h.Because of the implicitness of the problem implied by the dependence of Q nsol on T s , this correction is done iteratively during the computation of the TASFs.ERA-Interim SSV fields are used on a regular latitudelongitude grid with a resolution of 0.728 (500 3 251 points), and TASFs are calculated on the same grid.As illustrated in Fig. 1, q 2m is derived from the 2-m dewpoint temperature and the SLP [Eq.(B8); appendix B, section f].The saturation-specific humidity of air q sat is derived from the saturation partial pressure of water vapor e sat : q sat (T s , P) 5 «e sat (T s ) where « ' 0.62, P 0 is the SLP, and e sat (T s ) is calculated with the Goff-Gratch formula (Goff 1957) [Eq.(B1); appendix B, section a] as recommended by the World Meteorological Organization (WMO).We compute r from the equation of state of moist air, with T and q being adjusted from 2-to 10-m height within the bulk algorithm: where R d is the specific gas constant of dry air and r y ' 0.61.
We use instantaneous values of C P and L y , calculated according to Eqs. (B4) and (B5) (appendix B, sections c and d).TASFs are set to zero over ice-covered regions.For U z , we use the absolute 10-m wind speed vector.The effect of using U z relative to the sea surface current is tested in experiment CurrE (section 3e).
The 1982-2014 climatologies of wind stress and total turbulent heat flux Q T (i.e., Q L 1 Q H ) derived in CtrlE are shown in Fig. 2.
All the FORTRAN routines used to compute TASFs in this study are developed and maintained by the authors as part of the AeroBulk1 package and are currently being implemented as the default forcing function of the upcoming major release of the NEMO2 OGCM (version 4.0).
The approximations and assumptions studied in sections 3 and 4 are chosen to reflect many of the common parameterization issues faced by ocean and atmosphere practitioners when configuring bulk formulas for use in modeling experiments.As such, aspects of the bulk parameterization that are more esoteric in nature, and generally unknown to the practitioner, such as the choice of the stability functions or the Webb correction (Webb et al. 1980), are not considered in this study.They do, however, represent a source of uncertainty in computed TASFs and remain to be quantified in the manner presented here.

Impacts of approximations in the bulk formulas
We perform a series of sensitivity experiments to quantify the impact of various approximations in the bulk formulas [Eq.( 1)].For every sensitivity experiment, we FIG. 1. Schematic of the bulk approach used to compute turbulent air-sea fluxes based on the sea surface and near-surface atmospheric state variables.Note that T and q are not always provided at the same height as the wind speed (usually 2 and 10 m, respectively).Therefore, prior to using bulk formulas [Eq.( 1)], u and q are adjusted to the standard height of 10 m during the computation of transfer coefficients (in the bulk algorithm).
use the same configuration as in CtrlE (section 2b), except for the particular aspect being tested.The impact of the tested aspect is then assessed and discussed by analyzing global estimates of the differences in TASFs between the sensitivity experiment and CtrlE.

a. Density of air: Equation of state versus fixed value
In marine conditions, r [Eq.( 4)] typically spans values from 1.1 kg m 23 in the tropics to 1.3 kg m 23 in the cold environment, implying a variation of roughly 610% around the value of 1.2 kg m 23 often used as a reference3 (Table 2).Compared to the climatology of r computed in CtrlE using Eq. ( 4), a constant density of 1.2 kg m 23 overestimates r in the tropics by 3%-4% on average (Fig. 3).At higher latitudes, this underestimates r by about 2% in midlatitudes and 6% in subpolar regions.
In sensitivity experiment RhoE (Table 1), we use a fixed r value of 1.2 kg m 23 .An error in r affects all computed TASFs [Eq.( 1)].The spatial pattern of the mean relative error in wind stress and Q T (not shown) is almost identical to that of r (Fig. 3) but with the opposite sign.The momentum flux to the ocean increases by 3% in the tropics; however, because of relatively low winds, the resulting absolute effect is about 2 mN m 22 on average, and never exceeds 5 mN m 22 (Figs.4a,e).The error is small in subtropical highs where the year-round r is close to 1.2 kg m 23 .In mid-and high latitudes, the underestimation of r in RhoE leads to an underestimation of the wind stress with a mean difference being typically 210 mN m 22 (Fig. 4a) and the maximum differences (likely associated with winter storms) amounting to 230 mN m 22 (Fig. 4e).
In the tropics, RhoE overestimates E by 0.15 mm day 21 on average (equivalent to a heat flux deficit of 4 W m 22 ; Figs. 5a,e).In the extratropics, the increase in Q T , due to the underestimation of r during the cold season, dominates over less frequent decreases in Q T due to warm air episodes in summer.The resulting surplus of heat and freshwater input for the ocean at midlatitudes is about 13 W m 22 and 10.1 mm day 21 on average and is slightly more pronounced for Q T in the Northern Hemisphere (NH) due to the larger contribution of Q H . North of 408N, by limiting jQ T j in winter, RhoE leads to an underestimation of the amplitude of the seasonal cycle of Q T by almost 10 W m 22 (gray line in Fig. 5a).In coldair outbreak conditions (mostly dominated by Q H ), jQ T j can be underestimated by more than 20 W m 22 in RhoE (Fig. 5e).In February, the mean jQ H j over the Labrador Sea is underestimated by about 15 W m 22 (not shown).
The mean response of the global ocean to the fixed density assumption is dominated by the overestimation of evaporation in the tropics, leading to a cooling of 20.4 W m 22 and a freshwater deficit of 0.1 Sv (1 Sv 5 10 6 m 3 s 21 ; Table 1).

b. Accuracy of the estimate of q s
Inaccuracy in the estimate of q s has a direct effect on E and Q L , because of the linear dependence of E on the term Dq 5 q s 2 q z .It also has an indirect effect on other TASFs because of the dependence of all BTCs on the stability (which partly depends on Dq).Other indirect effects of an error in q s include slight modifications of the height-adjusted q z , r(z), and U B .Our results show that all these indirect effects are secondary and negligible compared to the direct error induced by Dq in Eq. (1c).

1) SALT-RELATED REDUCTION IN q s
Sud and Walker (1997) showed the importance of the salt-related 2% vapor pressure reduction in the estimate of q s in experiments with an AGCM.Zeng et al. (1998) and Brunke et al. (2002Brunke et al. ( , 2003)), using direct TASF measurements, concluded that omitting this correction leads to substantial errors in evaporation under strong wind conditions.
In sensitivity experiment SaltE (Table 1), the 2% reduction in Eq. ( 1f) is suppressed.The resulting increase in evaporation is about 5%-10% over the globe (Fig. 6a).These errors are quite outsized with respect to a correction of only 2%.The relative error in Dq, induced by the omission of the 2% correction, can be expressed (in %) as Err(Dq) 5 2 1 2 q z /q s .( 5) Equation ( 5) indicates that in marine conditions the relative error in Dq is always larger than 2% and that it has the potential to become substantial in the nearneutral to weakly stratified ASL (as q z /q s approaches 1).With T s 5 158C and RH 2m 5 85%, the relative error in Dq is 8% if T 2m 5 138C and 14% if T 2m 5 158C.
In SaltE, the regions that undergo the strongest increase in evaporation (larger than 10%) are strongly localized and are collocated with the regions with a TABLE 2. Impact of the (second row) temperature on the density and the (third row) specific heat capacity of moist air for a fixed SLP of 1010 hPa and a relative humidity of 80%.For the density, values in parentheses are for a SLP of 1040 (bold) and 970 hPa (italic).The relative errors are calculated as the difference between values obtained with a temperature of 308C or 258C and the reference value obtained with a temperature of 158C, divided by the same reference value.propensity for stably stratified ASLs, such as the eastern Pacific equatorial upwelling, the central midlatitude North Atlantic, and the Malvinas Current retroflection (Fig. 6a). Figure 7a shows that the error is the largest under weakly stratified conditions (DQ 5 118C), which leads to an overestimation of evaporation by 140%.Note, however, that under such conditions, E is small, so the large relative error may not necessarily imply large absolute errors.Therefore, it is mainly the regions with strong evaporation, and not necessarily those that exhibit large relative errors, that are affected the most, like the tropics (equatorial upwelling regions aside; Figs.5b,f).Nevertheless, the mean increase of jQ T j in SaltE is substantially larger than 2% over the whole global ocean (Fig. 6a).This highlights the outsized response of evaporation to a small change in q s [Eq.( 5)], particularly in marine conditions as q z is generally close to q s .Thus, for the global ocean, omitting the 2% reduction in q s leads to the relatively large excess of evaporation of 1.12 Sv, associated with a deficit of heat of 7.5 W m 22 (Table 1).The most likely error in Q T in the tropics is between 210 and 212 W m 22 , while it is twice as small (26 W m 22 ) in the extratropics of both hemispheres (Fig. 5f).The excess of evaporation is the highest at 6158 of latitude, where it leads to a deficit of freshwater and heat of 0.5 mm day 21 and 11 W m 22 , respectively, that decreases to zero toward the poles.
2) ACCURACY OF FORMULA FOR q sat The dependence of q sat on T s and P 0 is often approximated with a simplified formula instead of using the combination of Eqs. ( 3) and (B1).For instance, for the so-called CORE forcing used to produce ocean hindcasts and OGCM intercomparisons (Griffies et al. 2009), Large andYeager (2004, 2009) use the following formula: The relative error induced by the use of qsat instead of q sat from Eqs. ( 3) and (B1) never exceeds 62%, except in warm environments under stable conditions (Fig. 7b).
The error is minimized when r [in Eq. ( 6)] is calculated at the air-sea interface (at z 5 0, using T s and q s , suggested by the curve for T 2 T s 5 0 in Fig. 7b).If density is fixed, the relative error becomes substantially large outside of the 158-258C SST range (Fig. 7b).Therefore, the statement by Large and Yeager (2004, p. 7) that accuracy in q sat is ''not necessary given the uncertainties in the 0.98 factor'' is only consistent if Eq. ( 6) uses an accurate estimate of r [Eq.( 4); calculated at z 5 0 with T s and q s , preferably].However, most practitioners, when not using a fixed r, likely use Eq. ( 6) with r calculated at the height at which T and q are provided (2 or 10 m).Therefore, in sensitivity experiment QsatE (Table 1), we use the same approach: q s is calculated using Eq. ( 6) with r calculated with T and q at z 5 2 m [Eq.( 4)].
In QsatE, the largest errors in Q T are confined to the tropics, particularly over the equatorial upwelling region of the eastern Pacific where enhanced evaporation leads to a deficit of freshwater and heat flux for the ocean of 4% on average (not shown).This results in a zonally averaged heat flux deficit of about 4 W m 22 (Fig. 5c).In the extratropics, the underestimation of q s results in a small gain of heat for the ocean of less than 2 W m 22 on average (Fig. 5c).In the NH extratropics, the amplitude of the seasonal cycle of Q T decreases by about 5 W m 22 due to FIG. 3. Mean (1982-2014) relative error in the estimate of the air density due to the fixed density assumption (r 5 1.2 kg m 23 ), with respect to the reference density calculated with equation of state [Eq.( 4)].
the decrease and increase of evaporation in winter and summer, respectively.Globally, the enhanced evaporation in the tropics dominates over the extratropical evaporation reduction (Figs.5c,g), leading to a global surplus of evaporation of 0.18 Sv, equivalent to a heat flux deficit of 1.2 W m 22 for the ocean (Table 1).
The use of qsat together with the fixed density assumption should produce larger errors in evaporation (Fig. 7b).This is what our experiment QsatRhoE (Table 1) shows (Figs.5d,h); the errors are larger than those in RhoE and of opposite sign, due to reduced evaporation in the tropics (17 to 19 W m 22 on average) and enhanced evaporation at midlatitudes (roughly 25 W m 22 ).
The probability of errors exceeding approximately 5 W m 22 is also greatly increased due to the important contribution of the reduced evaporation in the tropics (cf.Fig. 5g and Fig. 5h).In QsatRhoE, the dominant contribution of reduced tropical evaporation yields a global deficit of evaporation of 0.21 Sv associated with a surplus of heat flux amounting to 1.7 W m 22 (Table 1).

c. SLP: Prescribed versus fixed value
The SLP is needed to compute r [Eq.( 4)] and q sat [Eq.( 3)].In many applications, such as ocean-only experiments, the use of a fixed SLP of 1010 hPa is common practice.Cyclone activity climatologies report the depth  1).Negative discrepancies mean a decrease of the heat flux to the ocean.White contours show the mean relative discrepancy in multiples of 5%.
of intense cyclone events being deeper than 960 hPa, with extreme cyclones approaching a central pressure of 920-930 hPa (Tilinina et al. 2013).The use of the fixed SLP assumption in such low pressure conditions has two effects: 1) q s is underestimated (Table 3), which implies a reduction of E and Q L with no significant impact on other TASFs (section 3b), and 2) r is overestimated, which enhances all TASFs.Under anticyclonic circulation conditions, the opposite scenario occurs.In the SLP range from 950 to 1040 hPa, both r and q s vary by roughly 65% (Table 3).
In sensitivity experiment SLPE, the SLP was fixed to 1010 hPa (Table 1).The fact that the rand the q s -induced errors have an opposite impact on E suggests a potential compensating effect.To quantify the individual contributions of these errors, an additional sensitivity experiment was run with the fixed SLP assumption implemented only in the calculation of r.
Over the Southern Ocean, strong cyclonic activity leads to an overestimation of r in SLPE, which enhances the wind stress by up to 5 mN m 22 over the Antarctic Circumpolar Current (ACC; Fig. 4b).North FIG. 7. (a) Scatterplot, in the U 10 , DQ space, of the relative error in Q T due to the omission of the salt-related 2% decrease in q s (experiment SaltE compared to CtrlE, a positive DQ means a stable ASL); only points for which the relative error is larger than 8% are shown, based on 6296 points in space that evenly spans the global ocean from 608S to 608N for years 2002-14 at a frequency of 3 h.(b) Relative error in q s , emerging from the use of the simpler formula qsat [Eq.( 6)] instead of Eqs. ( 3) and (B1), as a function of the SST, at a fixed SLP of 1010 hPa, for five different configurations of air-sea temperature difference. of 458N, periods of cyclonic activity alternate with blocking episodes more frequently, resulting in a smaller effect (12 mN m 22 ).No impact on the wind stress is found in the tropics and subtropics.
In SLPE, E and Q T are mainly affected by the q sinduced error in E (not shown), which is 3 times as large as the r-induced error on average and of opposite sign.The mean error in E and Q T does not exceed 60.1 mm day 21 and 62 W m 22 , respectively (Fig. 5i).Significant increases in Q T only occur over the ACC due to reduced evaporation (12.5 W m 22 is equivalent to 20.1 mm day 21 ).In the subtropics of both hemispheres, persistent anticyclonic conditions result in an overestimated q s , which enhances E and decreases Q T (20.1 mm day 21 and 22 W m 22 on average).In the NH, the decrease in Q T caused by high pressure is never less than 210 W m 22 , while in the Southern Hemisphere (SH), extreme low pressure episodes can increase Q T by 15 W m 22 (Fig. 5m).
In SLPE, the enhanced subtropical evaporation dominates over the reduced evaporation over the ACC.Consequently, the fixed SLP assumption enhances global evaporation by 0.1 Sv and leads to a global deficit of heat of 0.7 W m 22 (Table 1).
Implementation of the fixed SLP and fixed air density assumptions together, tested in experiment SLPRhoE, acts to increase the magnitude of the errors in E and Q T with respect to RhoE and SLPE (cf. Figs. 5a,i,j).In SLPRhoE, the fixed SLP assumption acts only through q s , and the compensating effect of the r-induced error that acts in SLPE vanishes.This amplifies the excess of evaporation in the tropics and subtropics and decreases Q T by 4-5 W m 22 between 308S and 308N (10.15 and 10.2 mm day 21 ).In the extratropics, the opposite effect occurs: an underestimation of E and an overestimation of Q T , particularly in subpolar regions (20.2 mm day 21 and 16 W m 22 at 608S).In northern mid-and subpolar latitudes, the substantial decrease in the amplitude of the seasonal cycle of Q T observed in RhoE is further decreased in SLPRhoE and reaches 210 W m 22 at 608N (cf.Figs.5a and 5j).

d. Cool skin and warm layer
Similar to COARE, the ECMWF algorithm is meant to be used with the SSST; it includes the same cool skin scheme as in COARE (Fairall et al. 1996a) and the warm layer scheme of Zeng and Beljaars (2005).In contrast to COARE and ECMWF, the NCAR algorithm is calibrated and developed for the (bulk) SST.
The warm layer scheme of the ECMWF algorithm is more advanced than its COARE counterpart, which tends to promote the diurnal warming of the surface layer induced by Q s , mostly under relatively calm conditions (Fig. 8a).Therefore, in regions with significant insolation and relatively weak winds, like the tropical Indian Ocean, the warm pool, and eastern boundaries, COARE tends to underestimate the mean effect of the diurnal warming of the surface layer (and therefore T s ) with respect to ECMWF, by about 0.18C (cf.Fig. 8b and Fig. 8c).Over the warm pool T s is 0.38C colder than the SST on average according to COARE, while it is 0.28C colder according to ECMWF.
In sensitivity experiment NoSkinE (Table 1), we suppress the COARE CSWL correction used in CtrlE, which implies that T s is set equal to the daily SST of ERA-Interim.In CtrlE, the SSST is colder than the SST by about 0.18 to 0.38C on average (Fig. 8b).In NoSkinE, the slightly warmer T s substantially enhances E and jQ T j (Fig. 6b); the resulting deficit of heat (and freshwater) for the ocean varies from roughly 10 W m 22 (0.3 mm day 21 ) in low latitudes to 2-3 W m 22 (0.03-0.06 mm day 21 ) in the subpolar regions (Fig. 5k).The contribution of Q H is not negligible and amounts to roughly 2 W m 22 from the equator to midlatitudes.The amplitude of the seasonal cycle of Q T is increased by up to 6 W m 22 in the NH extratropics (Fig. 5k) on account of calm and sunny conditions being more frequent than in the SH during summer.
Interestingly, the global distribution of DQ T is bimodal (Fig. 5o) due to two distinct responses to the skin correction between the tropics and extratropics.In the tropics, a low wind/high evaporation mode dominates, with a Gaussian distribution of the errors in Q T centered about 211 W m 22 .In midlatitudes, however, the magnitude of T s 2 SST is limited by the impact of the stronger winds (Fig. 8a) and also because the diurnal solar warming becomes more efficient in compensating the weaker jQ L j.
In summary, the omission of the skin correction in the COARE algorithm leads to an increase of global ocean TABLE 3. Sensitivity of the (second and third rows) sea surface saturation-specific humidity and (fourth and fifth rows) the density of moist air to the SLP.The density is calculated for a relative humidity of 80%.The relative errors are calculated as the difference between values obtained with a SLP of 1040 and 950 hPa divided by the reference value obtained with a SLP of 1010 hPa.evaporation by 0.74 Sv, associated with a deficit of heat for the ocean of 7 W m 22 (Table 1).

e. Wind vector referenced to the surface current
In many applications of bulk formulas, TASFs are determined from the apparent wind speed relative to the sea surface current vector u 0 , so that in Eq. ( 1), U z is replaced by U z 2 u 0 .
In sensitivity experiment CurrE (Table 1), we use the apparent 10-m wind vector of ERA-Interim relative to the prescribed daily u 0 from the Global Ocean Reanalysis and Simulation (GLORYS),4 mesoscale, eddypermitting ( 1 /48 resolution) ocean reanalysis (Ferry et al. FIG. 8. Temperature difference between the surface skin temperature T s and the bulk SST, calculated using the cool skin and warm layer schemes of the COARE and ECMWF algorithms.2012; Balmaseda et al. 2015).Note that CurrE is limited to the 1993-2013 period, which corresponds with the longest overlap in the GLORYS and ERA-Interim datasets.
The bulk formulas imply that the relative error in E, Q L , and Q H is of the order of u 0 /U, while it is at least twice as large for the wind stress on account of the quadratic dependence with wind speed.This suggests that the error is only significant if ocean currents are sufficiently strong with respect to the wind speed.This configuration predominates near the equator where CurrE shows significant relative errors in wind stress (Fig. 9a), particularly over the South Equatorial Current (SEC), the North Equatorial Countercurrent (NECC), and the North Equatorial Current (NEC).These results are consistent with previous studies reporting this effect to be the largest in the tropics (Luo et al. 2005;Dawe and Thompson 2006).
In CurrE, the wind stress decreases by about 2%-3% over the ACC (26 mN m 22 ), by 10% at the equator over the NEC (24 mN m 22 ), and by 2% in the NH extratropics (22 mN m 22 ; Fig. 4c).Globally, the decrease of wind stress ranges between 21 and 25 mN m 22 , with the strongest biases reaching 225 mN m 22 (Fig. 4g).The contribution from counterflow situations, in which the wind stress is increased, is significant and limits the overall decrease in wind stress by about 11%.Thus, the mean zonal error is entirely canceled over the NECC, on account of the ocean counterflow relative to the trade winds.
In CurrE, E and jQ T j are decreased in the tropics, which leads to a mean surplus of heat (and freshwater) flux of about 12 W m 22 (10.08 mm day 21 ) with a maximum over the SEC at the equator (15 W m 22 ; 11.9 mm day 21 ) and zero mean error over the narrow NECC band (Fig. 5l).The error in Q T diminishes outside of the tropics on account of the error being related to u 0 /U (Figs.5l,p), meaning that the surplus of heat flux is less than 1 W m 22 over the ACC as well as in the NH extratropics.Over the global ocean, errors in the surface heat flux very rarely exceed 5 W m 22 .The distribution is skewed toward positive error values (heat input to the ocean) but also displays significant contributions from the negative errors related to counterflow situations.
In summary, the correction of the wind speed with the daily u 0 leads to a global deficit of evaporation of 0.17 Sv and a global surplus of heat flux to the ocean of 1.2 W m 22 (Table 1).

f. C P and L y
The specific heat capacity of moist air C P [Eq.( 1b)] is slightly larger than its dry air counterpart C Pd (appendix B, section c).In typical marine conditions, C P is close to 1020 J kg 21 K 21 , with a range of 12.3%/21.2%[see Table 2 and Eq. ( B4)].Therefore, on a global scale, the errors in Q T due to the use of a fixed C P are negligible, given the small contribution of Q H to Q T (Fig. 2d).However, such an approximation could have a nonnegligible impact on episodes of extreme fluxes (Gulev and Belyaev 2012) associated with oceanic deep convection in subpolar and high latitudes, when most of the buoyancy loss to the atmosphere is due to Q H .In CtrlE, Q H can locally reach 2500 W m 22 in winter over the western branches of NH west boundary currents, the Labrador Sea, the Nordic Seas, and the southern ACC.
The latent heat of vaporization of water L y depends only on the temperature.Based on the widest spread of realistic marine conditions, the possible error in L y does not exceed 1.5% (Table 4).
Two separate experiments, using a fixed L y (2.46 10 6 Jkg 21 ) and C P (1005 J K 21 kg 21 ), respectively, confirm the negligible impact of these two approximations on the climatologies of Q T and E (not shown).The largest errors are observed in the fixed L y experiment, in particular in the tropics where enhanced evaporation leads to a mean decrease of Q T by 21.1 W m 22 (not shown).The Q T errors due to the fixed L y and C P approximations are of opposite sign in regions where they are the largest (tropics and extratropics, not shown).Consequently, experiment CpLvE, which combines both approximations, shows even smaller Q T errors (Table 1).
Therefore, both C P and L y can safely be assumed as constants.

g. Combined effects of common assumptions
We investigate to what extent multiple simplifications and approximations, implemented together, lead to accumulation or cancellation of errors in the time mean of computed TASFs.We address this by computing TASFs in experiment SimplE, which implements what we consider to be the most common choice of parameterization implemented in OGCMs.As such, SimplE differs from CtrlE in the following ways: d no CSWL correction of the SST (T s 5 SST), d assumed constant r (1.2 kg m 23 ), d estimation of q s using qsat [Eq.( 6)] rather than q sat [Eq.( 3)], and d assumed constant L y (2.46 J kg 21 ) and constant C p of dry air (1005 J kg 21 K 21 ).
Note that the second and third simplifications make the constant SLP assumption redundant.
Our results indicate that for the wind stress, only the impact of the fixed density assumption is expressed in SimplE (cf.Fig. 4a and Fig. 4d).
In SimplE, the errors in Q (Figs. 5q, 6c) are roughly equivalent to the addition of the errors from QsatRhoE FIG. 9. Mean relative wind stress discrepancy, with respect to the control experiment CtrlE, for experiments (a) CurrE, (b) C35E, (c) ECE, and (d) NCE (Table 1).Positive discrepancies mean an increase of the magnitude of the wind stress.White contours show the mean relative discrepancy in multiples of 5%.and NoSkinE, respectively (Figs. 5d,k).Thus, most of the tropical decrease in Q associated with the omission of the skin correction (NoSkinE) is compensated by the increase in Q due to the fixed density assumption used in conjunction with qsat (QsatRhoE).The net effect of these compensating errors leads to a relatively small error in Q in the tropics of approximately 22 W m 22 .In the extratropics, however, the approximations implemented in QsatRhoE and NoSkinE both induce enhanced evaporation and a decrease in Q.This results in a substantial enhancement of the extratropical heat and freshwater loss from the ocean in SimplE by about 6-8 W m 22 and 0.25-0.3mm day 21 , respectively.
The globally averaged errors in Q and evaporation in SimplE are 25 W m 22 and 10.54 Sv, respectively (Table 1).

Algorithm-related flux discrepancies
The bulk algorithm determines the BTCs and adjusts the temperature and humidity of air to the standard height (z 5 10 m).The algorithm may also include a convective gustiness contribution to the wind speed under calm and unstable conditions.Zeng et al. (1998) and Brunke et al. (2002Brunke et al. ( , 2003) ) carried out an extensive intercomparison of algorithms, including COARE 3.0 and the predecessors of ECMWF and NCAR algorithms, against direct flux measurements in different regions over the ocean.They reported a very large spread in the results and concluded that COARE and ECMWF were the least problematic of all.Since then, COARE, NCAR, and ECMWF have undergone updates and improvements, often taking advantage of direct flux measurements from recent observations.
To calculate the BTCs, every bulk algorithm relies on an empirical closure.In algorithms such as COARE and ECMWF, the roughness length z 0 is related to the friction velocity u* (Smith 1988): where n is the kinematic viscosity of air; a, the Charnock parameter (Charnock 1955), is typically 0.011 but varies substantially between studies [e.g., 0.035 in Garratt (1992)].Observational evidence has lead scientists to suggest that a increases with the wind speed in the range of moderate and strong winds (Yelland and Taylor 1996;Hare et al. 1999;Fairall et al. 2003).In COARE 3.0, a is set to 0.011 in the range from 0 to 10 m s 21 , linearly increases from 0.011 to 0.018 in the range 10 and 18 m s 21 , and remains equal to 0.018 for higher wind speeds (Fig. 10a).The ECMWF algorithm relies on a fixed a equal to 0.018 when IFS is not coupled to ECMWF's FIG. 10.Bulk empirical closures as a function of the wind speed at 10 m for COARE (red), NCAR (blue), and ECMWF (green) bulk algorithms.(a) Charnock parameter a; for NCAR, a is calculated according to Eq. (B10) (appendix B, section g).(b) Neutral drag and moisture transfer coefficients (thick and thin lines, respectively), as functions of the neutral wind speed at 10 m.  2013), aiming at improving the estimate of C D under strong wind conditions in the COARE algorithm, suggested a linear increase of a with the wind (U N10 ) up to a plateauing upper value of 0.028 for winds above 30 m s 21 (Fig. 10a).
The NCAR algorithm is based on a closure proposed by Large andPond (1981, 1982), in which the three neutral BTCs are empirical functions of U N10 .As noted by Beljaars (1997, p. 35), ''neutral transfer coefficients and surface roughness lengths are compatible concepts,'' so that the z 0 (u*) type of closure can be easily translated to the C N10 D (U N10 ) type [Eq.(B10); appendix B, section g].Similar to COARE 3.5, the implied equivalent a in the NCAR algorithm also suggests a roughly linear increase of a from moderate breeze to gale conditions, with a maximum value of 0.012 reached at a wind speed of about 25 m s 21 (Fig. 10a).The hurricane correction introduced in the latest update of the NCAR algorithm (Large and Yeager 2009) prevents the growth of C D under very strong winds, which implies a decreasing value of a for winds above 25 m s 21 .

a. Influence of the Charnock parameter
The value of a directly influences the growth of C D with the wind speed (Fig. 10).We assess the sensitivity of the computed wind stress to the value of a in two sensitivity experiments, C35E and ECE, in which C D is computed using COARE 3.5 and ECMWF algorithms, respectively.
Significant wind stress discrepancies between C35E and CtrlE are a result of the larger C D of COARE 3.5 under stormy conditions (typically U 10 .15 m s 21 ; Fig. 10b).Thus, in C35E, the wind stress is increased by about 15 mN m 22 over the ACC (Figs. 9b,4i) and 7.5 mN m 22 in the northern subpolar latitudes.Under extreme storm conditions, the local wind stress computed with COARE 3.5 can deviate from COARE 3.0 by more than 50 mN m 22 (Fig. 4m; not shown).The discrepancies between C35E and CtrlE are smaller in the North Atlantic midlatitudes on account of the lower probabilities of such high wind speeds.We also observe a slight decrease of the mean wind stress in the tropics due to the smaller C D in COARE 3.5 under calm and moderate winds.The associated deficit of momentum flux to the ocean amounts to typically 2 to 3 mN m 22 (Fig. 4i).
In ECE, the larger a under calm to moderate gale conditions (0.018 against 0.011 in CtrlE) leads to a globally uniform increase of the wind stress by 8%-9% on average (Fig. 9b).This results in an increase of the wind stress by 115 mN m 22 over the ACC, decreasing to about 3 mN m 22 at the equator (Fig. 4k).

b. Uncertainties related to the choice of the bulk algorithms
BTCs typically differ by about 10% between the bulk algorithms (Fig. 10).We note here that the COARE and ECMWF algorithms have similar behaviors, so we focus mainly on the discrepancies between the COARE and NCAR algorithms.We address this by performing sensitivity experiment NCE (Table 1), in which the NCAR algorithm (Large and Yeager 2009) is used to compute the BTCs and U B and adjust u and q from 2 to 10 m.Since the NCAR algorithm is intended to be used with a SST rather than a SSST (section 3d), TASFs are computed with T s 5 SST in NCE.
For wind speeds above 5 m s 21 , the C D of NCAR is smaller than that of COARE (Fig. 10b).This leads to a substantial reduction of the wind stress in NCE (Fig. 9c): from zero, on average, at the equator down to 230 and 217 mN m 22 over the ACC and northern midlatitudes, respectively (Fig. 4j).Under gale conditions and above (18-28 m s 21 ), the magnitude of the wind stress computed with COARE is typically 15% larger than that computed with NCAR (Fig. 10b).Therefore, the disagreement between NCE and CtrlE is the highest in the regions with the highest wind speeds, such as the ACC and the North Atlantic storm track (Fig. 9d).Under hurricane force winds, the latest version of the NCAR algorithm (Large and Yeager 2009) constrains C D to remain below 2.4 3 10 23 , while in COARE C D increases linearly with the wind speed (Fig. 10b).This can lead to large disagreements between the two algorithms under extreme wind conditions, which are illustrated by the heavy tail in the negative flux discrepancy range of the distribution (Fig. 4n).From calm up to light breeze conditions (U 10 , 5 m s 21 ), the C D of NCAR is larger than that of COARE (Fig. 10b).These conditions occur quite frequently (positive region of the distribution in Fig. 4n), but they have a limited and localized impact on the global wind stress disagreement in the time mean due to the small wind stress magnitudes involved, particularly outside of the tropical band (Fig. 9c).
Differences in evaporation and Q T emerge under all but gentle and moderate breeze conditions (3 , U 10 , 8 m s 21 , depending on the stability) on account of the larger C E (Fig. 10b) and generally larger C H (not shown) in NCAR compared to COARE (Fig. 10b).Under calm conditions (U 10 , 3 m s 21 ), C E and C H are also larger in NCAR than COARE.
The large disagreement in E and Q T between NCE and CtrlE (Figs. 5r, 6d) is the consequence of two effects: the enhanced E and jQ T j in NCE due to the larger C E and C H of NCAR and the decreased E and jQ T j in CtrlE due to the corrected (colder) T s .Thus, in NCE, Q T is decreased by about 12 W m 22 at the equator and 17 W m 22 in the subtropics (E is increased by 0.4 and 0.6 mm day 21 , respectively).In the extratropics, the contribution from Q H to this difference is not negligible compared to that from Q L and accounts for roughly onefourth of the difference (Fig. 5r).The mean relative differences in Q T are about 220% over the ACC and 215% in the northern midlatitudes (Fig. 6d).The absolute differences (not shown) are the largest over the west boundary currents (Gulf Stream and Kuroshio) where the magnitude of Q T is large (Fig. 2b).Over the Gulf Stream, daily discrepancies in Q T between NCE and CtrlE can locally amount to more than 300 W m 22 in winter.In NCE, the amplitude of the NH seasonal cycle of Q T strongly increases compared to CtrlE, with a maximum increase of 118 W m 22 at 408N (Fig. 5r).
The globally averaged differences in Q T and E between NCE and CtrlE are 213 W m 22 and 11.5 Sv.For reference, the total continental freshwater runoff into the ocean is estimated to be about 1.2 Sv (Lagerloef et al. 2010).When the COARE algorithm is used without the skin temperature correction, the agreement with NCAR is better, with the globally averaged differences between NCE and NoSkinE being 26.3 W m 22 and 10.8 Sv.

c. Treatment of calm conditions
In the unstable ASL, as the wind speed approaches zero, the production of TKE is dominated by buoyancy.Under such convective and calm conditions, standard bulk formulas are not suitable as they suggest that evaporation and turbulent heat fluxes tend to zero.These conditions were investigated in laboratory and model studies (Liu et al. 1979;Golitsyn and Grachev 1986;Godfrey and Beljaars 1991;Beljaars 1995), which led to the addition of a parameterized convective gustiness contribution to the wind speed in the bulk algorithm.The COARE and ECMWF algorithms both implement the same convective gustiness scheme.In the NCAR algorithm; however, the zero wind singularity is avoided by simply setting a minimum value for the scalar wind speed of 0.5 m s 21 .
To assess the impact of this simpler approach, we perform sensitivity experiment NoGustE (Table 1), similar to CtrlE, but in which very low wind speeds are dealt with using the same approach as in the NCAR algorithm.
The different treatment of extremely weak wind conditions has no significant impacts on the computed wind stress (Figs.4l,p).
In NoGustE, evaporation is slightly reduced over the global ocean (Figs.5t,x), especially in the tropics where the associated mean surplus of heat flux is about 1 W m 22 on average.The gustiness correction can significantly enhance U B with respect to jU z j in the range 0 to a few meters per second (not shown).As a consequence, E is decreased by about 4% over the warm pool and in Indonesia in NoGustE (not shown).Over the global ocean, the implementation of the simple minimum wind speed approach instead of the convective gustiness correction leads to a net warming of 11.1 W m 22 and a deficit of evaporation of 0.14 Sv (Table 1).

Summary and concluding remarks
We have quantified the impacts of some common simplifying assumptions and approximations in the parameterization of bulk formulas on climatologies of computed turbulent air-sea fluxes over a series of sensitivity experiments.Estimates of the 3-hourly wind stress, evaporation, and turbulent heat flux were calculated for the period 1982-2014 using the SST and near-surface atmospheric state variables from ERA-Interim.Each sensitivity experiment is evaluated with respect to our control setup, utilizing what we consider as the most consistent and appropriate parameterizations of bulk formulas.
The majority of approximations used in practice have limited impacts on computed flux errors.When tested individually, these errors only amount to about onequarter of the discrepancies related to the computation of the transfer coefficients themselves.However, when several approximations are implemented together, their summed contributions can lead to significantly larger errors, up to about 20% in the time and zonal mean.For example, the use of constant values for air density and sea level pressure (rather than those instantaneously computed) lead to flux errors of about 5% in subpolar regions when tested individually.These impacts combine when the approximations are implemented together, leading to a 12% increase of the heat flux over the ACC.
The accuracy of the estimate of the saturation-specific humidity (and consequently the temperature) at the airsea interface is found to be a major source of uncertainty in the computation of the surface turbulent heat flux and evaporation.The air above the sea is generally close to saturation, and, as such, a small error in the estimate of the sea surface saturation-specific humidity was demonstrated to have an outsized and amplifying effect on the evaporation errors.In particular, implementation of the salt-related 2% reduction of the saturation humidity and the consideration of the cool skin and warm layer are each shown to produce uncertainties of about 10 W m 22 and 0.3-0.4mm day 21 in the tropics, respectively.The use of a simplified formula for the saturation humidity is found to have weaker impacts, leading to errors no larger than 4%.However, large errors are possible when this simplified formula is used together with the constant air density assumption, leading to heat flux errors of 18 W m 22 in the tropics and 25 W m 22 in midlatitudes.
We note that bulk estimates of the wind stress are much less affected by approximations or uncertainties in state parameters than are the evaporation and turbulent heat fluxes.This is due to the more limited dependency of the wind stress on the stability of the atmospheric surface layer.
The largest uncertainties in the estimate of the turbulent air-sea fluxes are related to the computation of the transfer coefficient and therefore to the choice of the bulk algorithm.The comparison between the COARE 3.0 and NCAR algorithms, two widely used algorithms within the ocean community, shows large disagreements.NCAR is found to systematically enhance evaporation with respect to COARE by 12%-15% on average.This results in a mean extra cooling of 13 W m 22 and freshwater loss of 1.5 Sv over the global ocean.This cooling can reach up to 20 to 40 W m 22 over parts of the global ocean, especially in the ACC and in the large upwelling areas of the eastern basins.
These algorithm-related disagreements are found to be primarily linked to the empirical closure used to relate the transfer coefficients to the wind speed in neutral conditions.In particular, it is the value of the Charnock parameter, and more generally it is the value of the parameter that constrains the linear growth of the drag coefficient with the wind speed that matters.Disagreements related to these latter aspects amounts to about 10%, with the NCAR wind stress being notably weaker at latitudes higher than 308 in both hemispheres, especially over the ACC where it reaches 230 mN m 22 in the time and zonal mean.The latest improvements of the COARE algorithm (version 3.5) suggests that the drag coefficient of COARE 3.0 is likely too small in strong wind conditions, which in turn suggests an even greater tendency for NCAR to underestimate the wind stress, typically under gale force winds.
In the absence of surface gravity wave information from observations or wave model output data, such as spectral density, age, height, and steepness, the bulk, empirical, wind speed-dependent closures studied in this paper still represent the best choice for estimating the transfer coefficients.In that sense, the large disagreements between COARE and NCAR highlight the crucial need of further validation of these algorithms with respect to observations and new developments in their formulation.Therefore, efforts to converge on a more universally accepted dependence of the Charnock parameter or the neutral drag coefficient on the wind speed remain key priorities.

FIG. 2 .
FIG. 2. Mean climatology of (a) the wind stress module and (b) the total turbulent heat flux Q T (latent 1 sensible), obtained with the control experiment CtrlE.Corresponding zonal averages: (c) wind stress and (d) the total turbulent heat flux and its latent Q L and sensible Q H components and the net longwave radiative flux Q IR .

FIG. 4 .
FIG. 4. Wind stress errors, with respect to the control experiment CtrlE, for various sensitivity experiments (Table1).Errors are computed for the periods 1993-2014 in CurrE (denoted by *) and 1982-2014 for all other experiments.First and third rows indicate the time-mean and zonal-mean errors of the wind stress for all values (solid black lines), positive values only (red shading), and negative values only (blue shading).Positive values indicate a gain of momentum for the ocean.Second and fourth rows indicate the probability density functions of the 3-hourly time varying errors over the global ocean (608S-608N, gray shaded bars), tropics (238S-238N, red lines), and extratropics (608-308S in yellow and 308-608N in dashed blue, respectively).The distributions are computed using a bin size of 1 mN m 22 .
FIG. 8. Temperature difference between the surface skin temperature T s and the bulk SST, calculated using the cool skin and warm layer schemes of the COARE and ECMWF algorithms.(a) A function of the wind speed for two values of the surface downwelling solar radiative flux Q swY in the unstable (plain lines, DQ 5 22 K) and stable (dashed lines, DQ 5 12K) ASL; for a fixed relative humidity of 80% and Q lwY 5 350 W m 22 .(b) 1982-2014 mean with COARE.(c) 1982-2014 mean with ECMWF.
the formula of e sat (T ) given in Eq. (B1).g.Relation between a and C N10 D Charnock parameter a is expressed as a function of u* and z 0 from Eq. (7), and z 0 is deduced from the expression of the neutral drag coefficient at 10 m:

TABLE 4 .
Impact of sea surface temperature on the latent heat of vaporization of seawater.The relative errors are calculated as the difference between values obtained with a SST of 308 or 218C and the reference value taken at 158C, divided by the same reference value.